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We investigate a microscopic model for strongly correlated electrons with both on-site and near- 
est neighbor Coulomb repulsion on a 2D square lattice. This exhibits a state in which electrons 
undergo a "somersault" in their internal spin-space (spin-flux) as they traverse a closed loop in 
external coordinate space. When this spin-f/2 antiferromagnetic (AFM) insulator is doped, the 
ground state is a liquid of charged, bosonic meron- vortices, which for topological reasons are created 
in vortex-antivortex pairs. The magnetic exchange energy of the distorted AFM background leads to 
a logarithmic vortex-antivortex attraction which overcomes the direct Coulomb repulsion between 
holes localized on the vortex cores. This leads to the appearance of pre-formed charged pairs. We 
use the Configuration Interaction (CI) Method to study the quantum translational and rotational 
motion of various charged magnetic solitons and soliton pairs. The CI method systematically de- 
scribes fluctuation and quantum tunneling corrections to the Hartree-Fock Approximation (HFA). 
We find that the lowest energy charged meron- ant imeron pairs exhibit d-wave rotational symmetry, 
consistent with the symmetry of the cuprate superconducting order parameter. For a single hole in 
the 2D AFM plane, we find a precursor to spin-charge separation in which a conventional charged 
spin-polaron dissociates into a singly charged meron-antimeron pair. This model provides a uni- 
fied microscopic basis for (i) non-Fermi-liquid transport properties, (ii) d-wave preformed charged 
carrier pairs, (iii) mid-infrared optical absorption, (iv) destruction of AFM long range order with 
doping and other magnetic properties, and (v) certain aspects of angled resolved photo-emission 
spectroscopy (ARPES). 
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I. INTRODUCTION 

The microscopic understanding of the effect of charge 
carrier doping on spin— | antiferromagnetic ( AFM) Mott 
insulators is the central issue of the high-temperature 
superconducting cupratesJil Many puzzling experimental 
features of these systemsH suggest that a fundamental 
law of nature remains to be recognized. Extremely low 
doping ( S ~ 0.02 — 0.05 charge carriers per site) leads to 
a complete destruction of the long-range AFM order, and 
a transition to an unusual non-Fermi-liquid metal. This 
unusual metal becomes superconducting, with the transi- 
tion temperature T c strongly dependent on the doping S. 
The maximum T c is reached for dopings around S = 0.15. 
For higher dopings the critical temperature decreases to 
zero, and in the overdoped region a crossover towards 
a (non-superconducting) Fermi-liquid takes place. Two 
central questions require resolution. The first one con- 
cerns the nature of the charge carriers responsible for 
this non-Fermi-liquid metallic behavior. This is a fun- 
damental issue, since it lies outside the framework of 
Landau's Fermi-liquid theory and it necessitates under- 
standing the appearance of non-quasiparticlc-likc charge 
carriers in a system of interacting electrons. The second 
question concerns the nature of strong attractive pairing 
between these charge carriers, given the purely repulsive 
interaction between the constituent electrons. In conven- 
tional superconductors, the pairing attraction is due to 
overscreening of the electron-electron Coulomb repulsion 
by the ionic lattice. In the case of the high-tempexature 
superconducting cuprates, it has been suggested!!! that 
pairing is an intrinsic property of the electron gas itself 
mediated by AFM spin-fluctuations of the doped sys- 
tem. Accordingly, the challenge is to identify a strong 
attractive force based purely on repulsive Coulomb in- 
teractions. In this paper, we derive such a force and 
demonstrate that it leads to d-wave pairing of charge 
carrying holes introduced by doping a quantum, spin-i, 
Mott-Hubbard antiferromagnet. 

The simplest model Hamiltonians used to investigate 
the cuprate physics are the Hubbard model and the 
closely related t-J model. Unlike the ID problem, an 
exact solution for the 2D Hubbard Hamiltonian is not 
known. As a result, various approximations are neces- 
sary. Although the application of the mean-field the- 
ory has been severely criticized in this context, it pro- 
vides a valuable reference point for incorporating fluctu- 
ation effects. Moreover, even for the ID Hubbard model, 
essential features of the exact solution may be recap- 
tured by judiciously incorporating fluctuation and tun- 
neling effects into mean-field theory.0 The most straight- 
forward mean-field theory is the Hartree-Fock Approxi- 
mation (HFA). At half-filling (8 = 0) the HFA predicts 
an AFM Mott insulator ground-state. As the system 
is doped, HFA suggests that charge carrier holes in the 
AFM background assemble in chaxged stripes, which are 
quasi-one-dimensional structures .DO A large effort has 



been devoted to studying these charged stripes and re- 
lating them to certain features of the cuprates. 

Recently, a more fundamental investigation of the 
many-clcctron problem has suggested the possibility of an 
alternative model Hamiltonian for the cuprate physics. 
This model Hamiltonian, called the spin-flux model,El 
suggests that the long-range Coulomb interaction be- 
tween spin-i electrons leads to qualitative new physics, 
not apparent in the conventional Hubbard model ijsee 
Section 2). The results of the Hartree-Fock studyH of 
this spin-flux model are summarized in Section 3. They 
suggests that the undoped parent compound is also an 
AFM Mott insulator. However, unlike the conventional 
Hubbard model, the one-electron dispersion relations of 
the AFM mean-field in the spin-flux model match those 
measured experimentally through angle-resolved photo- 
emission for undoped cuprates. A proper description of 
the highest occupied electronic states (as provided by 
the spin-flux model), is crucial to considerations of dop- 
ing. The spin-flux model and the conventional Hubbard 
model differ dramatically in this regard. At the HF 
level, the doping holes added to the AFM background 
of the spin-flux model are trapped in the core of antifer- 
romagnetic spin vortices. This composite object (meron- 
vortex) is a bosonic charged collective mode of the many- 
electron system (the total spin of the magnetic vortex is 
zero). The reversal of the spin-charge connection pro- 
vides a microscopic basis for non-Fermi-liquid behavior. 

A magnetic vortex is strongly attracted to an antivor- 
tex. This attraction increases logarithmically with the 
distance between the vortex cores, and is stronger than 
the unscreened Coulomb repulsion between the charge 
meron- vortex cores. In effect, the increase in Coulomb 
energy between a given pair of holes is more than offset by 
the lowering in exchange energy between the background 
electrons as their vortices approach each other from far 
away. As the inter-vortex distance increases, more and 
more spins are rotated out of their AFM background ori- 
entation and the total energy of the system increases. 
Thus, even at the HF level, the spin-flux model provides 
a fundamental underpinning for the origin of both non- 
Fermi-liquid behavior, and strong pairing between the 
charge carriers. 

While providing a good starting point, the Hartree- 
Fock Approximation also has serious shortcomings. For 
instance, the ground-state wavefunction in the presence 
of doping is non-homogeneous (the static meron-vortices 
of the spin-flux model, or the charged stripes of the con- 
ventional Hubbard model, break translational symme- 
try). Physically, one expects that these charge carriers 
can move along the planes, resulting in a wavefunction 
which preserves the translational symmetry of the origi- 
nal Hamiltonian. Quantum dynamics of the charge car- 
riers also determines whether the doped ground-state is 
really a metal. Charge carriers in the optimally doped 
cuprates are quite mobile excitations, although their 
scattering rates are radically different from electrons in 
a conventional Fermi liquid. 
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A consistent way of treating the quantum dynamics of 
the charge carriers is ppesaded by the Configuration In- 
teraction (CI) Method,Ei3fl described in Section 4. Here, 
a linear combination of HF wavcfunctions is used in or- 
der to restore the various broken symmetries. For in- 
stance, in a doped system the CI wavefunction is chosen 
to be a linear combination of HF wavcfunctions with the 
charge carrier localized at different sites. Certain types 
of charge carriers can lower their total energy substan- 
tially by quantum mechanically hopping from one site 
to the next. We tested the accuracy of the CI method 
against the exact solutionis of the one-dimensional Hub- 
bard model in Reference 3. In the ID Hubbard model the 
CI method describes the quantum dynamics of charged 
domain wall solitons in the AFM background. By includ- 
ing these effects as fluctuation corrections to the Hartree- 
Fock mean-field theory, the CI method provides excel- 
lent agreement with the exact Bethe Ansatz solution for 
the ground-state energy of the doped ID Hubbard chain, 
over the entire U/t range. The CI method also leads 
to a clear demonstration of the spin-charge separation 
in ID. Addition of one doping hole to the half-filled an- 
tiferromagnetic chain results in the appearance of two 
different carriers: a charged bosonic domain- wall (which 
carries the charge but no spin) and a neutral spin- 1/2 
domain wall (which carries the spin but no charge). This 
studya demonstrates the effectiveness of the CI method. 
In this paper we use the CI method to investigate dynam- 
ics of the charged meron-vortices in the spin-flux model. 
Throughout this paper we exploit and refer to the anal- 
ogy between the charge excitations_of_the ID Hubbard 
model and the 2D spin-flux modelJI3't3 apparent in the 
CI approach.0 The CI results for the spin-flux model (pre- 
sented in Section 4) confirm that the meron-vortices are 
very mobile, suggesting that a collection of such mobile 
bosonic charge carriers is a non-Fermi-liquid metal. The 
CI method also allows us to identify the rotational sym- 
metry of the meron-antimeron pair wavefunction to be 
d-wave for the most stable pairs. An energetically more 
expensive metastable s-wave pairing is also possible. The 
possibility of spin-charge separation in 2D is elucidated. 
A summary of the results, their interpretation and con- 
clusions is provided in Section 5. 

II. THE SPIN-FLUX MODEL 

The effective 2D Hamiltonian that we use to describe 
the strongly correlated electrons residing in the 0(2p) — 
Cu(3d x 2_ y 2) orbitals of the isolated CuC-2 plane is the 
tight-binding model: 

n = -J2 (/'.• ( '- ( '- + h - c ) + V H n i n o C 1 ) 

id," i,j 

where ct creates an electron at site i with spin er, 
is the hopping amplitude from site j to site i on the 



square lattice, rii = ^2 c L c «<t is the total number of 

cr=l 

electrons at site i, and Vij is the Coulomb interaction 
between electrons at sites i and j. The dominant terms 
are the nearest-neighbor hopping tij = to and the on-site 
Coulomb repulsion Vu — U/2. If only these two terms 
are considered, and we shift the chemical potential by 
U, this reduces to the well-known Hubbard model. The 
neglect of all longer range Coulomb interaction (Vy = 0, 
if i ^ j), in the Hubbard model, is based on the Fermi- 
liquid theory notion of screening of the effective electron- 
electron interaction. However, Fermi liquid theory fails 
to explain many of the crucial features of the high-T c 
cuprates. In our description, we include the nearest- 
neighbor Coulomb repulsion, which we assume is on the 
energy scale of t. This leads to the generation of spin- 
flux, an entirely new type of broken symmetry in the 
many-electron system, which we show leads naturally to 
bosonic charge carriers in the form of meron-vortices, 
non- Fermi-liquid behavior and a strong attractive pairing 
force between holes in the AFM background. In order to 
extract the relevant physics from our starting Hamilto- 
nian, 

K = -hY, ( C ,W + h - c ) +Uj2n n n il +Vj2n i n j (2) 

we introduce bilinear combination of electron operators 
K 3 = c h a a0 c jP> M = M>2,3, for * ^ 3 (summation 
over multiple indexes is assumed). Here ct° is the 2x2 
identity matrix and a = (a 1 , a 2 , a 3 ) are the usual Pauli 
spin matrices. The notation means that the sites i 
and j are nearest neighbors. The quantum expectation 
value ( ) of the Ay operators are associated with charge- 
currents (n = 0) and spin-currents (/i = 1,2,3). Non- 
vanishing charge currents lead to appearance of electro- 
magnetic fields, which break the time-reversal symmetry 
of the Hamiltonian. Experimentally, this does not occur 
in the cuprates. In the following, we adopt the ansatz 
that there is no charge current in the ground state Ay- = 
but circulating spin-currents may arise and take the form 
Kj = ^Ay-n a , a = 1, 2, 3, where |Ay | = A for all i and 
j, and n is a unit vector. These spin-currents provide a 
transition state to the spin-flux mean-field that we use in 
this paper. 

Using the Pauli spin-matrix identity, |^(<''^/)* = 
5 aa i6ppi, it is possible to rewrite the nearest-neighbor 
electron-electron interaction terms as mnj = 2rii — 
^A^(A^) + . If we neglect fluctuations in the spin- 
currents, we can use the mean-field factorization, 
A&(A£)+- (A^)(A^)++A^<Ag)*-<Ag)(A5)*. Thus, 
the quartic nearest-neighbor Coulomb interaction term is 
reduced to a quadratic term that is added to the hopping 
term leading to the effective Hamiltonian: 

U = -* E (4aT%c j0 + h.c) "W- (3) 

<>,3'> i 

aft 
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Here, T^ g = (8 a p + iA^n ■ ff a p)/\fl + A 2 are spin- 
dependent SU{2) hopping matrix elements denned by 
the mean-field theory, and t — t yl + A 2 . In deriv- 
ing (|J) we have dropped constant terms which simply 
change the zero of energy as well as terms proportional 
to 53 n i which simply change the chemical potential. It 



was shown previouslyEHU that the ground state energy 
of the Hamiltonian of Eq.(|J) depends on the SU(2) ma- 
trices T lJ only through the plaquette matrix product 
T i2 T 23 T 34 T 4i = eX p(m • <r$). Here, $ is the spin- 
flux which passes through each plaquette and 2$ is the 
angle through which the internal coordinate system of 
the electron rotates as it encircles the plaquette. Since 
the electron spinor wavefunction is two- valued, there are 
only two possible choices for If $ = we can set 
T^ g — Sij and the Hamiltonian (^) describes conven- 
tional ordered magnetic states of the Hubbard model. 
The other possibility is that a spin-flux $ = tt pene- 
trates each plaquette, leading to T 12 T 23 T 34 T 41 = -1. 
This means that the one-electron wavefunctions are an- 
tisymmetric around each of the plaquettes, i.e. that as 
an electron encircles a plaquette, its wavefunction in the 
internal spin space of Euler angles rotates by 2tt in re- 
sponse to strong interactions with the other electrons. In 
effect, the electron performs an internal "somersault" as 
it traverses a closed path in the CuC>2 plane. This spin- 
flux phase is accompanied by a AFM local moment back- 
ground (with reduced magnitude relative to the AFM 
phase of the conventional Hubbard model) . In the spin- 
flux phase, the kinetic energy term in (^) exhibits broken 
symmetry as though a spin-orbit interaction has been 
added. However, it is distinct from the smaller, con- 
ventional spin-orbit effects which give rise to anisotropic 
corrections to superexchange interactions between local- 
ized spins in the AFM.E3 In the presence of charge car- 
riers this mean-field is unstable to the proliferation of 
topological fluctuations (magnetic solitons) which even- 
tually destroy AFM long range order. In this sense, 
the analysis which we present below goes beyond sim- 
ple mean field theory. The quantum dynamics of these 
magnetic solitons described by the Configuration Inter- 
action (CI) method, corresponds to tunneling effects not 
contained in the Hartree-Fock approximation. For sim- 
plicity, throughout this paper we assume that the mean- 
field spin-flux parameters T lJ are given by the simplest 
possible choice T 12 = -1,T 23 = T 34 = T 41 = 1 (see 
Fig. Q). In order to go beyond a mean-field description 
of the spin-flux, these matrices may also be treated as 
dynamical variables. In this paper, we go beyond mean- 
field theory in describing the antiferromagnetic degrees 
of freedom but restrict ourselves to a mean-field model 
of the spin-flux. 



III. HARTREE-FOCK RESULTS FOR THE 
SPIN-FLUX MODEL 

The Configuration Interaction Method utilizes a lin- 
ear combination- of judiciously chosen Hartree-Fock 
wavefunctions.Effl In this section, we provide a short re- 
view of the relevant Hartree-Fock results for the spin-flux 
model. A full comparison between the HFA for the spin- 
flux model and the conventional Hubbard model has been 
published elsewhere. □ 



A. The Static Hartree-Fock Approximation 

One of the most widely used approximations for the 
many-electron problem is the Static Hartree-Fock Ap- 
proximation (HFA). In this approximation the many- 
body problem is reduced to one-electron problems in 
which each electron moves in a self-consistent manner 
depending on the mean-field potential of the other elec- 
trons in the system. While this method is insufficient, 
by itself, to capture all of the physics of low dimensional 
electronic systems with strong correlations, it provides a 
valuable starting point from which essential fluctuation 
corrections can be included. In particular, we use the 
Hartree-Fock method to establish the electronic struc- 
ture and the static energies of various magnetic soliton 
structures. In the more general Configuration Interaction 
(CI) variational wavefunction, the solitons acquire quan- 
tum dynamics and describe large amplitude tunneling 
and fluctuation effects that go beyond mean field theory. 

In the HF approximation, the many-body wavefunc- 
tion |\&) is decomposed into a Slater determinant of ef- 
fective one-electron orbitals. The one-electron orbitals 
are found from the condition that the total energy of the 
system is minimized 



s lm =0 



(4) 



In order to approximate the ground state of the spin- 
flux Hamiltonian (|^), we consider a Slater determinant 
trial-wavefunction of the form 



I*) = 1141°)' 



(5) 



P =i 



where |0) is the vacuum state, N e is the total number of 
electrons in the system and the one-electron states are 
given by 



4 =X^™( i ' Cr ) C i<7 



(6) 



Here, the one-particle wave- functions <fi n (i,a) form a 
complete and orthonormal system. 

Using the wavefunction (|5|) in equation (^), and min- 
imizing with respect to the one-particle wavefunctions 



4 



<p n (i,a), we obtain the Hartree-Fock eigen-equations: 



E n cj) n (i,a) 



' £ J ■ 



U J2 (5^0(0 



0) (7) 



where (a x ,a y ,a z ) are the Pauli spin matrices and the 
charge density, 



iV, 



Q(») = (*l4 Cia |*) = y;i^,(i,a)i 



(8) 



and the spin density, 



S(i) 



(*l4^%l*> = ^ 



(9) 



must be computed self-consistently. The notation j S 
Vi appearing in (Q) means that the sum is performed 
over the sites j which are nearest-neighbors of the site i. 
The self-consistent Hartree-Fock equations (0,^,^) must 
be satisfied by the occupied orbitals p = 1 . . . N e , but can 
also be used to compute the empty (hole) orbitals. 

The ground-state energy of the system in the HFA is 
given by 



E, 



GS 



N ' ( 1 



S(i) 



(10) 



where the single particle energies are obtained from (Q). 

The approximation scheme described above is called 
the Unrestricted Hartree-Fock Approximation, because 
we did not impose constraints on the wavefunction 
which would require it to be an eigenfunction of various 
symmetry operations which commute with the Hamilto- 
nian (|^). If these symmetries are enforced, the method is 
called the Restricted Hartree-Fock Approximation. We 
use the Unrestricted HFA since it leads to lower ener- 
gies. The breaking of symmetries in our case implies 
that electronic, correlations are more effectively taken 
into account £3 The restoration of these symmetries is 
deferred until the CI wavefunction is introduced. 

In the undoped (half-filled) case, the self-consistent 
Hartree-Fock equations can be solved analytically for 
the infinite system, using plane-wave one-particle wave- 
functions. In the unrestricted Hartree-Fock approach, 
doping the system leads to the appearance of inhomoge- 
neous solutions, which break the translational invariance. 
In this case, we solve the unrestricted self-consistent 
Hartree-Fock equations numerically on a finite lattice. 
Starting with an initial spin and charge distribution S(i) 
and Q(i), we numerically solve the eigenproblem (ffl) and 



find the HF eigenenergies E n and wavefunctions (j) n (i, a). 
These are used in Eqs. (||) and @ to calculate the new 
spin and charge distribution, and the procedure is re- 
peated until self-consistency is reached. Numerically, we 
define self-consistency by the condition that the largest 
variation of any of the charge or spin components on any 
of the sites of the lattice is less that 10~ 9 between suc- 
cessive iterations. 



B. The undoped ground state 

For the undoped system, the Hartree-Fock equations 
(0) for an infinite system are easily solved. In the 
cuprates, long-range AFM order is experimentally ob- 
served. Accordingly, we choose a spin distribution at the 
site i — e x i x a + e y i y a of the form S(i) = Se, 
where e is the unit vector of some arbitrary direction, 
while the charge distribution is Q(i) = 1. In the spin- 
flux phase, it is convenient to choose a square unit cell, 
in order to simplify the description of the T lJ phase- 
factors. We make the simplest gauge choice compatible 
with the spin-flux condition for the T-matrices, namely 
that —T 12 = T 23 = T u — T 41 — 1 (see Fig. (§). This 
leads to the reduced square Brillouin zone —ir/2a < k x < 
7r/2a, —tt /2a < k y < ir/2a. From the Hartree-Fock equa- 
tions we find two electronic bands, characterized by the 
dispersion relations: 



E 



{ 8 f(k) = ±E sf (k) = ± 



2 f (k) + (US) 2 (11) 



where each level is four-fold degenerate and e s f(k) = 

—2t\J (cos (k x a)) 2 + (cos (k y a)) 2 are the nonintcracting 
electron dispersion relations in the presence of spin-flux. 
The HF ground-state energy of the spin-flux AFM back- 
ground is given by (see Eq. (|Io|)): 



E s J s = -AY J E sf (k) + N 2 u(s 2 + \ 



(12) 



where the AFM local moment amplitude is determined 
by the self-consistency condition M 



s=4£ us 



N 2 ^E sf {k) 



(13) 



At half-filling the valence band (E { sf } (fc) < 0) is com- 
pletely filled, the conduction band (E^'(k) > 0) is 
completely empty, and a Mott-Hubbard gap of magni- 
tude 2US opens between the valence and the conduc- 
tion bands. The ground-state of the undoped spin- flux 
model is an AFM Mott insulator. It is interesting to 
note that the quasi-particle dispersion relation obtained 
in the presence of the spin-flux (Eq. (|ll|)) closely resem- 
bles the dispersion as measured through angle-resolved 



5 



photo-emission spectroscopy (ARPES) in a compound 
such as Sr 2 Cu02Cl2cll (see Fig. |J). There is a a large 
peak centered at (tt/2, it/2) with an isotropic dispersion 
relation around it, observed on both the (0, 0) to (it, tt) 
and (0, tt) to (tt, 0) lines. The spin-flux model in HFA 
exhibits another smaller peak at (0, 7r/2) which is not 
resolvable in existing experimental data. This minor dis- 
crepancy may be due to next nearest neighbor hopping 
or other aspects of the electron-electron interaction which 
we have not yet been included in our modelJia The quasi- 
particle dispersion relation of the conventional Hubbard 
model (T 12 = T 23 = T 34 = T 41 = 1) has a large peak at 
(7r/2,7r/2) on the (0,0) to (tt, tt) line (see Fig. g), but it 
is perfectly flat on the (0, tt) to (tt, 0) line (which is part 
of the large nested Fermi surface of the conventional 2D 
Hubbard model). Also, it has a large crossing from the 
upper to the lower band-edge on the (0, 0) to (0, tt) line. 
This dispersion relation is very similar to that of the t — J 
model (see Ref. 17). 

While both the conventional and spin-flux model pre- 
dict AFM insulators at half-filling (at least at the HF 
level), the spin- flux also provides a much better agree- 
ment with the dispersion, relations, as measured by 
ARPES. As in the ID casej3 the effect of doping is the ap- 
pearance of discrete levels deep inside the Mott-Hubbard 
gap. These levels are drawn into the gap from the top 
(bottom) of the undoped valence (conduction) bands. 
Accordingly, the type of excitations created by doping 
depends strongly on topology of the electronic structure 
near the band edges. 

C. Charged solitons in the doped insulator: the 
spin-bag and the meron-vortex 

1. The spin-bag 

If we introduce just one hole in the plane, the self- 
consistent HFA solution is a conventional spin-polaron 
or "spin-bag" (see Fig. 3a). This tyge of excitation is 
the 2D analog of the ID spin-polaron a The doping hole 
is localized around a particular site, leading to the ap- 
pearance of a small ferromagnetic core around that site. 
The spin and charge distribution at the other sites are 
only slightly affected. In fact, the localization length of 
the charge depends on U jt, and becomes very large as 
US — > 0, since in this limit the Mott-Hubbard gap closes. 
For intermediate and large U jt, the doping hole is almost 
completely localized on the five sites of the ferromagnetic 
core. 

The spin-bag is a charged fermion, as can be seen by di- 
rect inspection of its charge and spin distributions. This 
is also confirmed by its electronic structure (see Fig. 3b). 
Thus, the 2D. spin-bag is indeed the analog of the ID 
spin-polaron.0 



2. The meron-vortex 

The 2D analog of the ID charged domain wall is the 
meron-vortex (see Fig. ^). Like the ID domain- wall, the 
meron-vortex is also a topological excitation, character- 
ized by a topological (winding) number ±1 (the spins 
on each sublattice rotate by ±27r on any closed contour 
surrounding the center of the meron). As such, a sin- 
gle meron-vortex cannot be created in an extended AFM 
background with cyclic boundary conditions by the intro- 
duction of a single hole (just as a single isolated charged 
domain wall cannot be created on an AFM (even) chain 
with cyclic boundary conditions, by the introduction of 
a single hole) . From a topological point of view, this is 
so because the AFM background has a winding number 
0, and the winding number must be conserved, unless 
topological excitations migrate over the boundary into 
the considered region. However, excitations can be cre- 
ated in pairs of total topological number 0. In the ID 
case, this means creation of pairs of domain walls, while 
in 2D this means the creation of vortex-antivortex pairs. 

From Figs. 4(a) and 4(b) we can see that the meron- 
vortex is a charged boson. The total spin of such a config- 
uration is zero, while it carries the doping charge trapped 
in the vortex core. Moreover, from its electronic spec- 
trum (Fig. 4b) we can see that only the extended states 
of the valence band are occupied. They are the only ones 
contributing to the total spin. Since only one state is 
drawn from the valence band into the gap, to become a 
discrete bound level, it appears that an odd (unpaired) 
number of states remains in the valence band. However, 
one must remember that for topological reasons, merons 
must appear in vortex-antivortex pairs. Therefore, the 
valence band has an even number of (paired) levels, and 
the total spin is zero. This argument for the bosonic 
character of the meron-vortex is identical tp-that for the 
charged domain wall in p^yacetyleneJijOEj 

Unlike in the ID case,0 we cannot directly compare 
the excitation energy of the spin-bag with the excita- 
tion energy of the meron-vortex. The reason is that the 
excitation energy of the latter increases logarithmically 
with the size of the sample, and therefore an isolated 
meron-vortex is always energetically more expensive than 
a spin-bag. However, topology requires that merons and 
antimerons are created in pairs. The excitation energy of 
such a meron- ant imeron pair is finite, allowing a mean- 
ingful comparison between excitation energies of a pair 
of spin-bags and a meron-antimeron pair. 

3. The meron-antimeron pair 

In Figs. 5(a) and 5(b) we show the self-consistent 
spin and charge distributions for the lowest energy self- 
consistent HF configuration found when we add 2 holes 
to the AFM background, in the spin-flux model, for 
U/t = 5. This configuration consists of a meron and an 
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antimeron centered on second nearest-neighbor sites. As 
a result of interactions, the cores of the vortices are some- 
what distorted. If the vortices were uncharged, vortex- 
antivortex pair annihilation would be possible. However, 
for charged vortices, the fermionic nature of the under- 
lying electrons prevents two holes from being localized 
at the same site, in spite of the bosonic character of the 
collective excitation. 

A very interesting feature of this configuration is the 
strong topological attraction between the vortex and the 
antivortex. The closer the two cores are to each other, 
the smaller is the region in which the spins are rotated 
out of their background AFM orientation by the vortices, 
and therefore the smaller is the excitation energy of the 
pair. Since the holes are localized in the cores of the vor- 
tices, this topological attraction between vortices is an 
effective attraction between holes in the purely repulsive 
2D electron system. This effect is unique to the spin- 
flux phase. In the conventional Hubbard model, vortices 
are not stable excitations. The vortex-antivortex attrac- 
tion increases as the logarithm of the distance between 
the cores. Therefore, the pair of vortices should remain 
bound even if full unscreened 1 jr Coulomb repulsion ex- 
ists between the charged cores, providing a compelling 
scenario for the existence of strongly bound pre-formed 
pairs in the underdoped regime. 

There is another possible self-consistent state for the 
system with two holes, consisting of two spin-bags far 
from each other (such that their localized wave functions 
do not overlap). The excitation energy of such a pair of 
spin-bags is simply twice the excitation energy of a sin- 
gle spin-bag. When this excitation energy is compared 
to the excitation energy of the tightly bound meron- 
antimeron pair, we find that it is higher by 0.15t (for 
U/t = 5). In fact, for U/t < 8 the HFA predicts that the 
meron-antimeron pair is the low-energy charged excita- 
tion, while for U/t > 8, the spin-bag is the low-energy 
charge carrier. This is analogous with the situation in 
ID, where the spin-bag was predicted to be^the low en- 
ergy excitation for U/t > 6.5, in the HFA.cl As in ID, 
however, we expect that this conclusion will be drasti- 
cally modified once the charged solitons are allowed to 
move along the planes and the lowering of kinetic energy 
through translations is also taken into consideration. 

We complete this review of the HF results by pointing 
out that the strong analogy between the ID Hubbard 
model and the 2D spin-flux model is due to the similar- 
ity between the electronic structures at zero doping. As 
seen from Fig.^, the 2D spin-flux model has isotropic dis- 
persion relations about the (7r/2,7r/2) point. This acts 
as a Fermi point for the noninteracting system as it does 
in the ID system. The two empty discrete levels drawn 
deep inside the Mott-Hubbard gap in the presence of the 
meron- vortex split from the (7r/2,7r/2) peaks of the elec- 
tron dispersion relation. The different topology of the 
large nested Fermi surface of the conventional Hubbard 
model leads to instability of the meron-antimeron pair. 
In fact, in the conventional Hubbard model doping holes 



assemble in charged stripes, as opposed to the liquid of 
meron-antimeron pairs which is the low-energy state of 
the doped spin-flux model. 



IV. CONFIGURATION INTERACTION 
METHOD RESULTS FOR THE 2D SYSTEM 

A. Configuration Interaction Method 

The essence of the Configuration Interaction (CI) 
method is that the ground-state wavefunction, for a sys- 
tem with N e electrons, is not represented by just a sin- 
gle N e x N e Slater determinant (as in the HFA), but 
a judiciously rGhosen linear combination of such Slater 
determinants. t2l Given the fact that the set of all pos- 
sible Slater determinants (with all possible occupation 
numbers) generated from a complete set of one-electron 
orbitals constitute a complete basis of the A^-particle 
Hilbert space, our aim is to pick out a subset of Slater 
determinants which captures the essential physics of the 
exact solution. 

Consider the CI ground-state wavefunction given by 



(14) 



where each \^>i) is a distinct N e x N e Slater determinant 
and the coefficients oli are chosen to satisfy the minimiza- 
tion principle: 
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i = l,N (15) 



This leads to the system of CI equations 

N N 



i=l J=l 



i=l,N (16) 



where E = (\P|7i |^ , )/( 1 J r |' l J r ) is the energy of the system in 
the \*f>) state , Hij = i\Ti.\^f j) are the matrix elements 
of the Hamiltonian in the basis of Slater determinants 
= 1, N}, and (9y = arc the overlap ma- 

trix elements of the Slater determinants (which are not 
necessarily orthogonal). The CI solution is easily found 
by solving the linear system of equations ( |l6[) , once the 
basis of Slater determinants {| , F;),z = 1,N} is chosen. 

If we denote by <f)p(i,cr) the p = 1, ...,N e one-electron 
occupied orbitals of the Slater determinant |\t n ), these 
matrix elements are given by: 
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The matrix elements of the Hamiltonian (Q) can be writ- 
ten as: 
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where the expectation values of the hopping and on-site 
interaction terms arc: 
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We now consider the specific choice of the Slater de- 
terminant basis (l^j),? = 1,N}. Strictly speaking, one 
may choose an optimized basis of Slater determinants 
from the general variational principle: 



o 



l,N;p=l,N e (19) 



However, implementation of this full trial-function min- 
imization scheme (also known._as a multi-reference self- 
consistent mean-field approachlij) is numerically cumber- 
some even for medium-sized systems. Instead, we select 
the Slater determinant basis = 1,N} from the 

set of broken symmetry, Unrestricted Hartree-Fock wave- 
functions (||), their symmetry related partners and their 
excitations. Clearly, (||) satisfies ([l9|) by itself, provided 
that the a coefficients corresponding to the other Slater 
determinants in Eq.(|l4|) are set to zero (see Eq. (f|)). 
Since this Unrestricted HF wavefunction is not trans- 
lationally invariant (the doping hole is always localized 
somewhere on the lattice), we can restore the transla- 
tional invariance of the CI ground-state wavefunction by 
also including in the basis of Slater determinants all the 
possible lattice translations of this Unrestricted HF wave- 
function. Furthermore, if the self-consistent configura- 
tion is not rotationally-invariant (e.g. a meron-antimeron 



pair), all possible rotations must be performed as well. 
By rotation we mean changing the relative position of 
the meron and antimeron while keeping their center of 
mass fixed. 

Clearly, all the translated HF Slater determinants lead 
to the same HF ground-state energy ( i f> n \7-l\ i & n ) = Eqs 
as defined by Eq. (10). The CI method lifts the degen- 
eracy between states with the hole-induced configuration 
localized at different sites (see Eq. (pf)), thereby restor- 
ing translational invariance. We may identify the lower- 
ing in the total energy due to the lifting of this degen- 
eracy as quantum mechanical kinetic energy of decon- 
finement which the doping-induccd configuration saves 
through hopping along the lattice. In addition, quan- 
tum fluctuations in the internal structure of a magnetic 
soliton can be incorporated by including the lowest or- 
der excited state configurations of the static Hartree- 
Fock energy spectrum. Such wavefunctions are given by 
a p ah\^}, where p > N e labels an excited particle state 
and h < N e labels the hole which is left behind (see 
Eq. (||)). Once again, all possible translations (and non- 
trivial rotations) of this "excited" configurations must 
be included in the full CI wavefunction. These additions 
can describe changes in the "shape" of the soliton as it 
undergoes quantum mechanical motion along the plane. 

The CI method is described in more detail in Refer- 
ence 3, where it is used to study the ID Hubbard chain 
in order to gauge its accuracy by comparing its results 
with the exact Bethe ansatz solution. We showed that 
the CI method recaptures the essential physical features 
of the exact solution of the ID Hubbard chain, such as 
spin-charge separation, as well as leading to remarkable 
agreement of ground state energies of doped chains for 
all values of U/t. The main difference between the ID 
case and the 2D case is the computation time required. 
The computation time for one matrix element Tl nm scales 
roughly like ./V 9 , where N is the number of sites. The 
number of configurations included in the CI set scales as 
N\/N S \(N — N s )l when N s solitons are present. For both 
an iV-site chain and an N x N lattice, the HF "bulk" 
limit is reached for N > 10. In the ID casea we used 
chains with N = 10 — 25, and numerical calculations 
can be easily performed. However, in 2D the smallest 
acceptable system has 100 sites, leading to an enormous 
increase in the computation time. Nevertheless, our sam- 
ple of results in 2D suggest a simple and clear physical 
picture which we describe below. 



B. Spin-bag Dissociation: 
Spin-Charge Separation in 2D 



The charged spin-bag carries a spin of 1/2. Let 
\^+), \^-) be the HF determinants for the spin-bag cen- 
tered at any two nearest neighbor sites, respectively, and 
let S z = J2i S z (i) = I J2i,<j crc L c «>- be tne total spin op- 
erator in the z-direction. Then, S z \^+) = while 
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s z \*. 



— i | \I/ ) (or viceversa), since moving the cen- 
ter of the spin-bag by one site leads to a flip of its total 
spin (see Fig. 3(a)). Consequently, (\I r _|\f f + ) = 0. Since 
the Hubbard Hamiltonian commutes with S z , it follows 
that = 0. From the CI equation (|l|) we 

conclude that there is no mixing between states with dif- 
ferent total spin. As a result, it is enough to include in 
the CI set only those configurations with the spin bag 
localized on the same magnetic sublattice. Let us denote 



by |* 



(0,0), 



the initial static Hartree-Fock configuration, 



and by |^( nm )) the configuration obtained through its 
translation by n sites in the x-direction and m sites in 
the y-direction (cyclic boundary conditions are imposed). 
The condition that only configurations on the same sub- 
lattice are included means that n + m must be an even 
number, and the cyclic boundary conditions mean that 
0<n<N-l,0<m<N-l, for a, N x N lat- 
tice. As explained in detail in the ID analysis!] mixing 
configurations with the charged spin-bag localized at dif- 
ferent sites and then subtracting out the contribution of 
the undoped AFM background allows us to calculate the 
dispersion band of the spin-bag itself: 

E sb (k) = E(k 1 N)-N 2 e GS . (20) 

Here, the total energy of the lattice with the spin-bag 



and the CI wave-function 



|%) = ^2 exp ( i ( k ^ n + k v\ m ) a )^{n,m)) 
(n,m) 

are the solutions of the CI equations (fl6|). The finite size 
of the lattice and cyclic boundary conditions restricts the 



calculation to fc-points of the form fc 
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Na 



ae x 



0e v ), 



where (a, 0) is any pair of integer numbers. As usual, 
only fc-points inside the first Brillouin zone need to be 
considered. 

An analysis of the dependence of the spin-bag disper- 
sion relation E s b(k) on the size N x N of the lattice is 
shown in Fig. |g[ for the conventional Hubbard model 
(upper panel) and spin-flux model (lower panel), and 
U/t = 5. We used 6x6, 8x8, 10x10 and (only for the 
spin-flux model) 12x12 lattices. The dispersion relation is 
plotted along lines of high symmetry of the full Brillouin 
zone. For comparison, we also show the excitation energy 
E^ b F obtained in the static HFA as a full line. For both 
models, we see that the spin-bag dispersion band is al- 
most converged, even though we used quite small lattices. 
The convergence is somewhat slower in the spin-flux case, 
as seen most clearly at the (0,0) point. Although the val- 
ues obtained from the four lattices all differ at (0,0), the 
extremum values correspond to the 6x6 and the 8x8 lat- 
tices, while the values for the 10x10 and 12x12 lattices 
are indistinguishable. We conclude that the fit ( p0| ) is 
legitimate. 



From Fig. || we also see that the dispersion relations 
for the spin-bag in the two different models are very dif- 
ferent. The dispersion relations over the full 2D Bril- 
louin zone are shown in Fig. 0, and they are seen to 
mimic the electronic dispersion relation of the underly- 
ing undoped AFM background, shown in Fig. [|. This is 
consistent with the quasi-particle nature of this charged 
spin- 1/2 spin-bag. In the conventional Hubbard model, 
the undoped AFM background has a large nested Fermi 
surface along the (0, n) to (n, 0) line, and it is exactly 
along this line that the spin-bag dispersion band has a 
minimum. Similarly, the lowest energy of the spin-bag of 
the spin-flux model is at (tt/2, tt/2), corresponding to the 
Fermi points of the underlying undoped spin-flux AFM 
background. 

The extra kinetic energy B a j,(^, §) — E^ b F saved by 
the spin-bag through quantum hopping is 0.37t in the 
conventional model and 0.56i in the spin- flux model (for 
U/t = 5). Since the spin-bag is confined to one mag- 
netic sublattice, it must tunnel two lattice constants to 
the next allowed site. Consequently, the energy gained 
through hopping (of order t 2 /U) is small. This is dis- 
played, for the spin-bag of the spin-flux model, in Fig. 
|8|, where we plot the lowering in kinetic energy of the 
deconfincd spin-bag E s b(%,^) — E^ b F as a function of 
t/U . A similar dependence for the spin-bag of the con- 
ventional Hubbard model is presented elsewhere.^?] As in 
the ID case, we conclude that the spin-bag in 2D is a 
rather immobile quasiparticle-like excitation. 

In the ID model it is energetically favorable for the 
immobile spin-bag to decay into a charged bosonic do- 
main wall and a neutral fermionic domain wall, resulting 
in spin-charge separation.^! The analog of the ID charged 
bosonic domain wall is the 2D charged bosonic meron- 
vortex of the spin-flux model. If the spin-bag decays into 
a charged meron- vortex, a magnetic antivortex must also 
be created for topological reasons. Unlike the pair of 
domain walls in the ID case, the vortex-antivortex pair 
is tightly bound by a topological binding potential that 
increases as the logarithm of the vortex-antivortex sep- 
aration. Therefore, we expect that the doping charge 
is shared between the two magnetic vortices. One tech- 
nical problem for testing this hypothesis is that such a 
configuration ( a vortex-antivortex pair sharing one dop- 
ing hole) is not self-consistent at the static Hartree-Fock 
level. In the static approximation we require two doping 
holes to stabilize two vortex cores and create a meron- 
antimeron pair. We can, however, construct a trial wave- 
function to describe the singly-charged vortex-antivortex 
pair, by adding one electron in the first empty state of the 
self-consistent doubly-charged meron- antimeron configu- 
ration. The first empty levels of the meron-antimeron 
pair are the localized levels bound in the vortex cores, 
two for each vortex (see Fig. 4(b)). Because of degener- 
acy between the two lower discrete levels of the pair, we 
have in fact two distinct trial wave- functions, obtained by 
adding one electron in either of these two lower localized 
gap electronic states of a self-consistent meron-antimeron 







pair. These wavefunctions are not invariant to rotations 
( see Fig. ||) . Therefore we must include in the CI set of 
Slater determinants the configurations obtained through 
7r/2 rotations of the vortex-antivortex pair about its fixed 
center of mass in addition to translated configurations. 
As a result, we have a total of 8iV 2 configurations describ- 
ing the singly charged vortex-antivortex pair localized at 
all possible sites with all possible orientations about the 
center of mass. 

We performed this CI analysis for a 10 x 10 lattice and 
U/t = 5. The HF energy of a simple static spin-bag is 
— 0.82i (measured with respect to the HF energy of the 
undoped AFM background, equal to — 76.76t). The en- 
ergy of the static singly-charged vortex-antivortex pair 
is — 0.23<. Thus, we see that because this singly-charged 
pair trial wavefunction is not self-consistent, in the static 
case this configuration is energetically much more costly 
than the self-consistent spin-bag configuration. However, 
if we allow for quantum motion of these configurations, 
the situation changes dramatically. Performing the CI 
analysis for the set of all possible translated spin-bag 
configurations, we find that the energy of the spin-bag 
is lowered to —1.2At. Performing the CI analysis for the 
set of all translated and rotated singly-charged vortex- 
antivortex pairs we find that this configuration's energy is 
lowered to — 2.18i. This shows that the vortex-antivortex 
pair has lowered its translational and rotational kinetic 
energy by almost 2t, thereby becoming the low-energy 
charge carrier. This large number is not surprising, since 
unlike the spin-bag, the vortex-antivortex pair is not con- 
strained to motion on one magnetic sublattice. As a re- 
sult, such configurations lower their kinetic energy by an 
amount on the scale of t, as opposed to t 2 /U for the 
spin-bag configuration. For larger U/t values this effect 
is even more pronounced. 

We conclude that these results strongly support the hy- 
pothesis of spin-bag dissociation into a much more mobile 
singly-charged vortex-antivortex pair, analogous to the 
ID spin-bag dissociation into a pair of a charged bosonic 
domain-wall and a neutral fermionic domain-wall. u Un- 
like in the ID case, however, we do not have distinct 
charge and spin carriers for the composite excitation. In- 
stead the spin and the charge are shared equally between 
the vortex and the antivortex. If, on the other hand, 
there was a mechanism whereby the vortices became 
unbound, complete spin-charge separation could occur, 
in which one vortex traps the hole (and is therefore a 
charged me ro jj i) .and the other vortex carries the spin, in a 
lotus-flowciil3H3 (or undoped magnetic meron) configura- 
tions. In the absence of the corresponding self-consistent 
static HF configuration we are not able to settle this ques- 
tion. At very low doping, the strong vortex-antivortex 
topological attraction binds the spin and charge together. 
This is different from the ID case, where the absence of 
long-range interactions between the domain-walls allow 
for a complete spin-charge separation at any doping and 
even at zero temperature. 

This scenario opens a new avenue for research into how 



the system evolves with doping. If each hole is dressed 
into a singly-charged vortex-antivortex pair, when two 
such pairs overlap it is possible that both doping charges 
move to the same pair, creating a meron- ant imeron pair 
of charged bosons. Such pre-formed charge pairs may 
condense into a superconducting state at low tempera- 
tures. The other uncharged vortex-antivortex pair may 
either collapse and disappear (this is likely to happen 
at low-temperatures) or remain as a magnetic excitation 
of the system (at higher temperatures), mediating the 
destruction of the long-range AFM order, the renormal- 
ization of the spin-wave spectrum and the opening of the 
spin pseudogap. 

C. D-wave pairing of charge carriers 

From the static HF analysis we found that the most 
stable static self-consistent configuration with two dop- 
ing holes added to the AFM background of the spin-flux 
model is the meron-antimeron pair, for 3 < U/t < 8. At 
larger U/t, two charged spin-bags become more stable, 
in the static HF approximation. This is in close anal- 
ogy to the prediction that the spin-bag is energetically 
more favorable than the static charged domam-\ssall for 
U/t > 6.5, in the HFA of the ID Hubbard model! How- 
ever, in the ID case the charged domain-wall is consid- 
erably more mobile than the charged spin bag, gaining a 
kinetic energy of the order of t as opposed to t 2 /U energy 
gained by the spin-bag. As a result, when this kinetic 
energy of deconfinement is taken into account within the 
CI method, the charged domain-wall is found to be the 
relevant charged excitation for all values of U/t. A sim- 
ilar picture emerges in the 2D case, because the meron- 
vortices are much more mobile than the spin-bags. 

For the 2D system, we have shown that the charged 
spin-bag has very similar behavior to the ID charged 
spin-bag. The analog of the ID charged bosonic domain- 
wall is the 2D charged bosonic meron-vortex. We now 
consider the properties of the doubly-charged meron- 
antimeron pair. All the numerical results quoted in the 
rest of this section refer to a meron-antimeron pair on a 
10x10 lattice, in the spin-flux model with U/t = 5. 

As already discussed, the meron-antimeron pair is not 
rotationally invariant. We can find the rotational ki- 
netic energy saved by the pair as it rotates about its 
center of mass. In the present case, only 4 configurations 
need to be included, corresponding to the four possible 
self-consistent arrangements of the meron and antimeron 
about their fixed center of mass (see Fig. ||). Simple rota- 
tion by 7r/2 of the one-particle orbitals <j> p (i, a) about the 
center of mass is not, however, sufficient to generate the 
rotated configurations. First of all, the ir/2 rotation also 
changes the spin-flux parameterization. If the spin-flux of 
the initial configuration is T 12 = -1, T 23 = T 34 = T 41 = 
1, a 7r/2 rotation leads to a state corresponding to the 
rotated configuration T 12 = 1,T 23 = -1,T 34 = T 41 = 1. 
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Thus, following the ir/2 rotation, a unitary transforma- 
tion must be performed in order to restore the initial 
spin-flux parameterization. For the case cited above, 
this simply implies the change in the one-particle orbitals 
4> p {i Xl i v , u) — > —4>p(ix, iyi f) for &U sites (i Xl i y ) which are 
a type '2' site of the unit cell, in other words sites with 
even i x and odd i y (also, see Fig. |l|). The second ob- 
servation is that the rotation by 7r/2 also changes (flips) 
all the spins of the AFM background surrounding the 
pair. Thus, an extra n rotation about an axis perpendic- 
ular to the lattice plane is necessary to restore the align- 
ment of the AFM background. Following these trans- 
formations it is straightforward to generate the Slater 
determinants l^), l^) and 1^4) corresponding to the 
meron-antimeron pairs rotated by ? , n and ^ from the 
initial self-consistent HF meron-antimeron pair described 
by The CI method can be used to find the rota- 

tional energy saved by superposing these rotated meron- 
antimeron configurations. The lowest CI energy found 
is 0.46£ below the energy of the static pair, and corre- 
sponds to d-wave symmetry. By this we mean that the 
coefficients multiplying the 4 rotated states in the CI 
wave-function |\&) = J2i=i a i\^i) satisfy the condition 

Oil = — Qf2 = Q!3 = —Ui. 

Translation of a pair over the whole lattice can also 
be investigated. Since the pair does not carry any spin, 
all possible translations must be included (there is no 
restriction to same magnetic sublattice configurations). 
This leads to a total of iV 2 possible configurations for 
a N x N lattice. Again, when various configurations 
are generated from the initial self-consistent HF meron- 
antimeron state care must be taken to preserve 
the same spin-flux parameterization and the same AFM 
background orientation. This can be achieved perform- 
ing transformations similar to the ones described above. 
As a result of performing the CI method on the set of 
translated states, we find the dispersion relation of the 
(un-rotated) meron-antimeron pair. This is shown in Fig. 
[u| In this plot we show the total energy of the lattice 
with the moving meron-antimeron pair, as a function of 
the total momentum of the pair. Quantum hopping of 
the meron-antimeron pair lowers its total energy by an 
extra 1.29t. Two other interesting features are observed 
in Fig. The first one is that the dispersion relation 
of this rigidly polarized pair is not invariant to rotations 
by 7r/2, as expected. More important is that the min- 
ima of the dispersion relation occur at the (jr, n) points. 
Since the momentum of the pair is twice the momentum 
of either the meron or the antimeron, this is consistent 
with the fact that in their lowest energy state, both the 
meron and the antimeron have momenta of (7r/2,7r/2), 
in the spin-flux model. The doubling of the size of the 
Brillouin zone is also a direct consequence of this dou- 
bling of total momentum ( for comparison with undoped 
dispersion relation, see Fig. ||). 

However, to obtain the true energy of the charged 
pair, we must mix all AN 2 rotated and translated meron- 
antimeron configurations. All have the same static HF 



energy and are equally important in the CI method. Let 
I ^0(0, 0)) denote the initial self-consistent static Hartree- 
Fock meron-antimeron configuration, and \tyg(n,m)) de- 
note the configuration obtained through translation by n 
sites in the x-direction and m sites in the y-direction, as 
well as a rotation by an angle of 6% of the pair about its 
center of mass. Here < 9 < 3 and < n < jV — 1, < 
m < N — 1 (cyclic boundary conditions are imposed). 
The CI wavefunctions are then given by: 



l*iE> 



^ ^ ^i(k x n-\-k y 7n)a 



> e |* 9 (n,m)) (21) 



(n,m) 



The dispersion relation 



-Epazr (&) 



E 



HF 



obtained from this complete set is shown in Fig. O. The 
reference point is the HF energy of the self-consistent 
meron-antimeron pair E%£ T — — 78.52i. Thus, we see 
that the total kinetic energy saved by the freely moving 
and rotating meron-antimeron pair is 1.75t, when the to- 
tal momentum of the pair is (tt,tt). This energy equals 
the sum 0.46t + 1.29i of rotational and translational ki- 
netic energies found before (the number of significant fig- 
ures indicates the estimated accuracy of the computa- 
tional method). The rotational invariance of the disper- 
sion band is also restored. Besides the absolute minima 
about the (tt, tt) points, there is a more shallow minimum 
region about the (0,0) point. 

The rotational symmetry of the meron-antimeron pair 
wavefunction, defined by the coefficients (ot$), is a func- 
tion of the total momentum carried by the pair, as shown 
in Fig. [l2[ The absolute minima points (tt, tt) and the 
area around them correspond to pairs with d-wave sym- 
metry. By this, we mean that the coefficients ag have 
the form 



ag = cxp 



K) 



n: 



(22) 



with J = 2, i.e. ao = —oti = cti = — Q3. The core area, 
about the local minimum (0, 0) point, corresponds to s- 
wave symmetry. In this region the coefficients (ag) again 
satisfy Eq. (p2[), but for J — 0, i.e. ao = a± — a-i = 
The intermediary area appears to be a mixture of differ- 
ent J values. A simple decomposition of the form ( [22] ) 
is no longer possible. Instead, a sum of such terms cor- 
responding to different J values is required. Since we 
only have rotations by tt/2, a unique identification of the 
composite symmetry is not possible. Moreover, the en- 
ergy of the states in this intermediary area is at the top 
of the dispersion band. In order to find the correct CI 
states for energies well above the static HF value (i.e. 
larger than zero, in this case) we must add to the CI set 
the first set of excited HF states. For a meron-antimeron 
pair, excitation of an electron from the valence band onto 
the empty localized levels inside the Mott-Hubbard gap 
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costs about 1.5t of energy, for U/t — 5, so such states 
should contribute significantly in the CI states with pos- 
itive energies and modify their dispersion and symmetry 
(for this reason, we do not show the upper three high- 
energy bands in Fig. ^) . Consequently, both the energy 
and the symmetry of the states in the intermediary area 
may be modified from the ones shown in Fig. How- 
ever, the minima at (n, ir) and (0, 0) are at energies well 
below zero. Their energies and rotational symmetry are 
unaffected by additions of higher energy configurations 
to the CI set. 

The fact that we obtain two distinct minima is not 
very surprising. As argued before, we expect that in- 
dividual merons and antimerons are created with mo- 
menta of (±7r/2, ±7r/2). As a result, two different cou- 
plings are possible. A (7r/2,7r/2) meron can pair with 
a (7r/2,7r/2) antimeron, creating a pair of total momen- 
tum (-7T, 7r). This is the most stable coupling, leading 
to the lowest possible energy of — 1.75£ below the static 
HF energy. This pair has d-wave symmetry The sec- 
ond possible coupling is between a (tt/2, tt/2) meron and 
a (— 7r/2,— 7r/2) antimeron. This pair has a total mo- 
mentum of (0,0), and s-wave symmetry. However, this 
coupling is less strong. For the U/t = 5 case considered, 
the energy of the s-wave (0,0) pair is 1.28i above the 
energy of the d-wave (n, it) pair. The existence of both 
d-wave and s-wave pairing, and the dominance of the d- 
wave pairing, have beep established experimentally for 
the high-T c cuprates.Eil We are not aware of any other 
microscopic theory that predicts the two types of pairing 
to appear in different regions of the Brillouin zone. 

The total kinetic energy saved by the meron-antimeron 
pair through quantum hopping and rotation is of order t, 
as expected, since the pair is not restricted to one mag- 
netic sublattice and tunneling is not required for motion. 
Consequently, we expect that the energy saved by the 
meron-antimeron pair for larger values of U/t is com- 
parably large. On the other hand, the energy saved by 
the spin bag through tunneling motion scales like t 2 /U. 
In fact, we argued that a spin-bag may dissociate into 
a singly-charged vortex-antivortex pair in order to en- 
hance its mobility. However, even if dissociation does 
not occur, the kinetic energy saved by a pair of spin-bags 
is significantly smaller than the kinetic energy saved by 
the meron-antimeron pair. This shows that for U/t = 5 
the meron-antimeron pair is even more favorable ener- 
getically than the HFA predicts and suggests that the 
U/t range where meron-antimeron pair formation occurs 
may extend well beyond the U/t = 8 limit found within 
the HFA.EI In the ID case, the range of stability of the 
charged domain wall versus the charged spin-bag is ex- 
tended, (from the HF prediction of U/t = 6.5) to all U/t 
range .0 A numerical analysis is needed to determine if 
the limit is extended to infinity in the 2D case as well. 



V. SUMMARY AND CONCLUSIONS 

The Configuration Interaction Approximation incorpo- 
rates crucial quantum translational and rotational mo- 
tion of the charge carriers, which are absent in the static 
Hartree-Fock Approximation. Given the accuracy of the 
CI method in recapturing certain features of the exact 
Bethe ansatz solution of the ID Hubbard model,B we be- 
lieve that the CI method is likewise a very powerful tool 
for describing effects beyond mean-field theory in 2D. For 
2D systems, numerical calculations are much more time 
consuming. However, our small sample of results is quite 
suggestive of a simple physical picture. In direct analogy 
with the ID results, we find that the bosonic charged 
meron- vortices are much more mobile than the fermionic 
charged spin-bags. The extra kinetic energy gained by 
the meron- vortices very likely extends their region of sta=, 
bility beyond the U/t = 8 limit suggested in the HFA.EI 
There are also strong indications that a charged spin- 
bag is unstable to decay into a singly-charged vortex- 
antivortex pair, analogous to the spin-charge separation 
in the ID case. Nucleation of such pairs of vortices with 
doping is expected to further influence the magnetic be- 
havior of the cuprates. The bound state and the unbound 
continuum states of the singly charged meron-antimeron 
pair may account for the anomalous "quasiparticle" spec- 
tral widths observed on angle-resolved photo-emission ex- 
periments. 

The symmetry of the meron-antimeron pairs emerges 
very clearly from the CI treatment. We find two regions 
of stability of the meron-antimeron pair. Pairs with to- 
tal momentum of (n, tt) have d-wave symmetry, and are 
the most stable. Pairs with total momentum (0, 0) have 
s-wave symmetry and have a smaller gap. Thus, we find 
that different pairing appears in different regions of the 
Brillouin zone. These results appear to have a direct 
bearing on numerous experiments, which show a mixture 
of strong d-wave component and a smaller s-wave com- 
ponent in the superconducting state of the cuprates E3 

Many other features of our model are in agreement 
with experimentally observed properties of the cuprate 
superconductors. Nucleation of magnetic vortices with 
doping explains a variety of magnetic properties, starting 
with complete destruction of the long-range AFM order 
for very low doping concentration. As we can see from 
Fig. 5(a) a tightly-bound meron-antimeron pairs disturbs 
the long-range AFM ordering of most of the spins on the 
10x10 lattice. For very low dopings, these pairs are far 
from each other, and there are many spins on the plane 
whose orientations are not affected by any pair. Thus, 
most of the spins maintain the long range AFM order. 
However, as the doping increases to about 5% (which is 
roughly equivalent to having two meron-antimeron pairs 
on the 10x10 lattice) the areas occupied by each meron- 
antimeron pairs start to overlap with those occupied by 
the neighboring pairs. At this doping the orientation of 
all the spins on the Cu02 planes is affected by at least 
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one pair of vortices, and therefore the LRO is severely dis- 
rupted. The local ordering, however, is still AFM. This 
picture explains the extremely low doping necessary for 
the disappearance of LR AFM order, as well as the fact 
that the spin correlation length is basically equal to the 
average distance between holes (vortices) and does not 
depend strongly on the temperature. E3 Each hole carries 
its vortex with it, and the spins in each vortex are cor- 
related with each other. The correlation length is thus 
roughly equal to the average size of the vortex which is 
defined by the average inter-vortex (inter-hole) distance. 
The nucleation of magnetic vortices also explains the split 
of the (tt, it) AFM Bragg peak into the four incommen- 
surate peaks whose positions shift with doping,ELas ob- 
served in LaCuO and, more recently, in YBaCuOo The 
form factor of any given vortex already gives rise to an 
apparent splitting of the neutron scattering peak. As 
demonstrated in Ref. 9, even at the mean-field level we 
recapture the neutron scattering data using a random 
distribution of meron-vortices. 

Optical behavior of the cuprates is also explained nat- 
urally using our model. Two features develop in the 
optical absorption spectra with doping: a broad mid- 
infrared temperature-independent absorption band, and 
a stpangly temperature-dependent low-frequency Drude 
tail.c3 In our model the broad mid-infrared band is re- 
lated to excitation of electrons from the valence band 
onto the empty levels bound in the vortex coresjj which 
are localized deep inside the Mott-Hubbard gap (see Fig. 
4(b)). The number of localized levels scales with the 
number of vortices, and inter-vortex interactions lead to 
their splitting into a broad band. As such, this mech- 
anism is similar to the one leading to a broad mick 
infrared absorption band in polyacetylene with doping.E£l 
The polyacetylene band is due to electronic excita- 
tions inside the cores of the polyacetylene domain-wall 
solitons, — rwhich are the topological analogues of meron- 
vortices.113113 Another strong argument in favor of this 
interpretation, is provided by photoinduced absorption 
experimentsEH If the undoped parent compounds are il- 
luminated with intense visible light, they develop absorp- 
tion bands that resemble the mid-infrared bands of the 
doped compounds. Similar behavior is observed in poly- 
acetylene, and is attributed to the nucleation of solitons 
by photoexcited electron- hole pairs .C3 The second com- 
ponent of the optical spectrum is the Drude tail. It 
results from the response of the freely moving charged 
vortices to the external electric field. The strong tem- 
perature dependence of this tail is determined by the 
scattering mechanism for merons (presumably due to in- 
teractions with other merons and spin- waves). This in- 
terpretation is also supported by the fact that the super- 
conducting transition leaves the mid-infrared absorption 
band unchanged. Merons with internal electronic struc- 
ture are still present on the planes but pair condensation 
leads to a collapse of the Drude tail into a 8(u>) response. 

As already discussed, nucleation of charged meron- 
vortices with doping provides a microscopic basis for 



non-Fermi-liquid behavior, due to the bosonic nature of 
the mobile charged excitations. The model also predicts 
the existence of pre-formed pairs with d-wave symmetry, 
which ate thought to be responsible for the pseudo-gap 
effects.c£l As the number of pairs of charged bosons in- 
creases with doping and the temperature is lowered, the 
meron-antimeron pairs Bose condense and become co- 
herent, leading to superconductivity. This mechanism 
naturally explains the puzzling scaling of the superfluid 
density with doping p s ~ 5, in other words with the num- 
ber of holes, not of electrons. In our model, it is the 
doping-induced positively charged meron-vortices that 
are the mobile charge carriers. As a result, the den- 
sity of preformed meron-antimeron Cooper pairs is ob- 
viously proportional to doping. Finally, for large dop- 
ings (S > 0.30 — 0.40) the average vortex size become 
extremely small and the very cores of the merons start 
to overlap. In this limit the Mott-Hubbard gap is com- 
pletely filled in by the discrete levels, and the spin-flux 
state becomes energetically unstable relative to a normal 
Fermi liquid. □ 

It is noteworthy that all of the independent features 
described above are in qualitative agreement with our 
model, which has essentially no free or adjustable pa- 
rameters. The choice of U/t is fixed by the experimen- 
tally measured size of the Mott-Hubbard charge transfer 
gap at zero doping. More detailed comparisons with spe- 
cific experiments may require the incorporation of spe- 
cific (smaller energy scale) interactions which are not in- 
cluded in this simplest version of the spin-flux Hamilto- 
nian. A derivation of the explicit consequences of this 
picture appears to be worthy. A more comprehensive 
and quantitative e comparison with the experiments may 
prove quite fruitful. 
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FIG. 1. Choice of the gauge for describing the mean- field 
spin-flux background. Physical observables depend on the ro- 
tation matrices T 1J only through the plaquette matrix prod- 
uct t 12 T 23 T 34 T 41 . Shown above is the simplest (spin inde- 
pendent) gauge choice describing a 2-7r-rotation of the internal 
coordinate system of the electron (described by 3 Euler an- 
gles) as it encircles an elementary plaquette. This is a new 
form of spontaneous symmetry breaking for a strongly inter- 
acting electron system in which the mean-field Hamiltonian 
acquires a term with the symmetry of a spin-orbit interaction. 
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FIG. 2. A comparison between the experimentally de- 
termined E(k) quasi-particle dispersion relation, from an- 
gle resolved photoemission studies (ARPES), for the insu- 
lating Sr2Cu02Ci2 (open circles with error bars) and the HF 
AFM spin-flux model dispersion relation (full line) and the 
HF AFM conventional Hubbard model dispersion relation 
(dashed-dotted line). Three directions in k— space are shown: 
(0,0) to (7T,7r), (tt,0) to (0,0) and (tt, 0) to (0,tt). While the 
peak on the (0,0) to (vr,7r) is equally well described in both 
models, the mean-field spin-flux model gives a much better 
agreement for the (tt,0) to (0,0) and (tt,0) to (0, tt) direc- 
tions. The fitting corresponds to U = 2.01 eV, t = 0.29 eV 
for the spin-flux phase, and U =1.98 eV, t = 0.21 eV in the 
conventional phase. The experimental results are the ARPES 
results of Ref. 17. 
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FIG. 3. (a) Self-consistent spin distribution of a 10x10 
lattice with a spin bag centered at (5,5). The spin-bag has 
a small ferromagnetic core, and the magnetic order is only 
locally affected, (b) Electronic spectrum of a spin bag on a 
10x10 lattice, for U/t = 5 in the spin-flux model. Eigenen- 
ergies E a are plotted as a function of a = 1, 200(= 2N 2 ). 
Only the first N 2 — 1 = 99 states are occupied. There are two 
empty bound discrete levels deep into the Mott-Hubbard gap 
(a — 100, 101), one of which comes from the valence band of 
the undoped AFM compound (see inset). There is also an 
occupied discrete level below the valence band (a = 1). The 
valence band is spin paired, since it has an even number of 
levels. Thus, the total spin of the spin-bag comes from the 
discrete occupied level. The spin-bag is a charged, spin-i 
fermion. 
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FIG. 4. (a). Self-consistent spin distribution of a 10x10 
lattice with a meron-vortex in the spin-flux phase. The core 
of the meron is localized in the center of a plaquette, in the 
spin-flux phase (in the conventional phase, the core of the 
meron-vortex is localized at a site). This excitation has a 
topological winding number 1, since the spins on either sub- 
lattice rotate by 2n on any curve surrounding the core. The 
magnitude of the staggered magnetic moments is slightly di- 
minished near the vortex core but is equal to that of the 
undoped AFM background far from the core. (b). Electronic 
spectrum of a meron-vortex on a 10x10 lattice, for U/t = 5, in 
the presence of the spin-flux. Eigenenergies E a are plotted as 
a function of a = 1, 200(= 2N 2 ). Only the first iV 2 - 1 = 99 
states are occupied (the valence band). There are two dis- 
crete empty levels deep into the Mott-Hubbard gap, one of 
which (a = 100) comes from the valence band of the undoped 
AFM parent. Merons must be created in vortex-antivortex 
pairs (for topological reasons). Each pair removes two lev- 
els from the undoped AFM valence band. Thus, the valence 
band remains spin paired, and the total spin of this excitation 
is zero. This meron is a spinlcss, charged, bosonic collective 
excitation of the doped antiferromagnet. 




FIG. 5. (a). Self-consistent spin distribution for a tightly 
bound meron-antimeron pair. The meron (M) and the an- 
timeron (A) are localized on neighboring sites. The total 
winding number of the pair is zero. The magnetic AFM order 
is disturbed only on the region where the vortices are local- 
ized. The attraction between holes is of topological nature 
and on long length scale is stronger than unscreened Coulomb 
repulsion between charges, (b). Self-consistent charge distri- 
bution for a tightly bound meron-antimeron pair. The dop- 
ing charge is mostly localized on the two plaquettes contain- 
ing the meron and antimeron cores. The two holes local- 
ized in the vortex cores are responsible for the fact that the 
meron-antimeron pair does not collapse (due to Fermi statis- 
tics, it is impossible to have two holes at the same site). 
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FIG. 6. Dispersion relation of the charged spin-bag E s b(k) 
(in units of t) plotted along lines of high symmetry in the 
Brillouin zone. The upper plot shows the dispersion band 
of the spin-bag in the conventional model, while the lower 
one shows the dispersion band of the spin-bag in the spin-flux 
model. U/t — 5 in both cases. Circles, squares, diamonds and 
triangles show the results obtained from CI analysis of 6x6, 
8x8, 10x10 and 12x12 lattices. We conclude that the results 
are already almost converged, even for such small lattices. 
The full lines show the excitation energy of the spin-bag at 
the static HF level. 
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FIG. 7. Dispersion band of the spin-bag in the conven- 
tional model (upper panel) and the spin-flux model (lower 
panel), for U/t — 5. We show the full 2D Brillouin zone (k is 
measured in units of 7r/a). The spin-bag dispersion relations 
have the symmetry as the dispersion relations of the under- 
lying undoped AFM background, shown in Fig. |^. While the 
(n/2a, n /2a) point corresponds to the minimum excitation en- 
ergy of the spin-flux spin-bag, in the conventional model all 
points along the (0,7r) to (7r,0) have almost the same energy. 
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FIG. 8. The extra kinetic energy (in units of t) gained by 
the spin-bag AE — E sb (^, ^) — Ef b F (circles) as a function of 
t/U. The log- log graph clearly shows the linear dependence. 
This is expected, since the spin-bag must tunnel two sites to 
the next allowed position. This is a second order hopping pro- 
cess and therefore this charged excitation is rather immobile. 







FIG. 9. Schematic representation of the four different 
meron-antimeron pair configurations obtained through n/2 
rotation about their fixed center of mass. The upper-left 
picture is a schematic representation of the self-consistent 
meron-antimeron pair shown in Fig.|H[ 




FIG. 10. Total energy E(k) (in units of t) vs. k of the 
lattice with an (un-rotated) meron-antimeron pair as a func- 
tion of the total momentum k of the pair. The momentum 
units are n/a. The HF energy of the static meron-antimeron 
pair is — 78.52i. Quantum hopping lowers the overall energy 
of the pair by 1.29t. Since the meron-antimeron configuration 
is not rotationally invariant, the dispersion relation is also not 
invariant to tt/2 rotations. 
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FIG. 11. The lowest energy dispersion band E pa i r (k) (in 
units of t) as a function of the total momentum k of the 
meron-antimeron pair. The momentum units are n/a. For 
convenience, the reference energy is taken to be the static HF 
energy of the self-consistent meron-antimeron pair. Quantum 
hopping and rotation lowers the overall energy of the pair by 
1.75t. The rotational symmetry of the dispersion relation is 
restored (compared to Fig. |Io| ). 
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FIG. 12. The rotational symme- 

try of the meron-antimeron wave-function as a function of 
the total momentum carried by the pair (measured in units 
of n/a). The outside region (containing the absolute minima 
points (7r,7r)) has d-wave symmetry (J = 2), while the core 
region about the (0,0) point has s-wave symmetry (J = 0). 
The intermediary area is a mix of s+d wave symmetry. 
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